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ABSTRACT 

An  efficient  modeling  of  the  electromagnetic  scattering 
properties  of  symmetric  lamellar  periodic  structures  at 
normal  incidence  is  discussed  within  the  framework  of 
the  modal-field  approach.  By  taking  advantage  of  the 
inherent  symmetry,  a  numerical  computational  speedup 
of  up  to  8  times  can  be  achieved  by  explicitly  ignoring 
anti-symmetric  modes,  which  play  no  role  in  the  scat¬ 
tering  process  at  all. 

1  INTRODUCTION 

We  consider  here  the  special  problem  of  electro¬ 
magnetic  scattering  of  normal  incident  waves  onto  a 
lamellar  periodic  structure  whose  unit  cell  possesses 
a  plane  about  which  the  structure  exhibits  a  reflec¬ 
tion  symmetry.  This  problem  is  commonly  encoun¬ 
tered  in  corrugated  quantum  well  infrared  photodetec¬ 
tors  (C-QWIPs).  (Choi  1997,  Mao  et  al.  2002,  Choi 
et  al.  2003,  2004a, b,c)  Other  applications  include  anti- 
reflection  coatings  and  interference  filters. 

QWIPs  have  applications  in  many  different  areas, 
including  industrial,  medical,  geological,  environmen¬ 
tal,  homeland  security,  and  defense.  In  the  case  of  de¬ 
fense,  they  are  critical  to  all  phases  of  any  ballistic  mis¬ 
sile  system,  including  missile  seeking  and  detection,  and 
decoy  countermeasure.  They  also  have  applications  in 
long  range  surveillance  and  land-mine  detection 

An  example  of  a  symmetric  C-QWIP  structure  is 
shown  in  Fig.  1.  The  structure  has  a  semi-infinite  sub¬ 
strate  below  and  a  semi-infinite  cover  region  at  the 
top,  with  3  grating  layers  sandwiched  between  them. 
These  grating  layers  consist  of  a  multiple  quantum- 
well  (MQW)  layer  of  thickness  tq  sandwiched  between 
two  semi-conductor  layers  of  thicknesses  tc  and  tb  .  The 
MQW  region  is  assumed  to  have  a  biaxial  dielectric 
constant  given  by  ev  in  the  vertical  direction  and  by  ex 
in  the  horizontal  direction.  The  width  for  the  air  region 


is  s,  and  that  for  the  MQW  region  is  ic,  and  so  the  lat¬ 
tice  spacing  a  is  s  +  w.  The  incident  wave  is  assumed 
to  be  a  plane  TM  polarized  wave  incident  from  below 
and  traveling  perpendicularly  to  the  layers. 
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Figure  1:  An  example  of  a  C-QWIP  structure. 


Modeling  and  optimizing  the  performance  of  plro- 
todetectors  like  the  present  one  presents  certain  chal¬ 
lenges.  First  of  all  there  must  be  strong  diffraction 
of  the  waves  within  the  structure.  This  means  that  the 
various  dimensions  of  the  structure  must  be  comparable 
to  the  wavelength  of  the  incident  radiation,  and  the  di¬ 
electric  contrast  must  also  be  sufficiently  strong.  Mod¬ 
eling  of  diffraction  in  such  a  strongly  scattering  medium 
is  known  to  be  always  challenging  and  requires  accurate 
handling  of  the  entire  set  of  Maxwell’s  equations  with 
its  associated  set  of  boundary  conditions. 

Secondly,  photodetector  structures  are  intrinsically 
multilayered.  The  above  simple  example  of  a  C-QWIP 
consists  of  only  5  layers  (3  periodic  layers  and  2  semi¬ 
infinite  homogeneous  layers).  More  realistic  devices  of¬ 
ten  have  metallic  layers  and  other  dielectric  layers  for 
various  electrical  and  mechanical  and  structural  pur¬ 
poses.  The  number  of  layers  quickly  increases  when 
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multi-color  operations  are  required.  The  handling  of 
periodic  structures  that  have  slanted  sidewalls  often  re¬ 
quires  the  artificial  sectioning  of  such  a  layer  into  multi¬ 
ple  sub-layers,  which  further  increase  the  total  number 
of  layers  that  must  be  handled  in  modeling  such  struc¬ 
tures.  It  is  also  well-known  that  many  methods  for 
imposing  the  required  electromagnetic  boundary  con¬ 
ditions  suffer  increasing  degree  of  instabilities  with  in¬ 
creasing  number  of  layers. 

Thirdly,  realistic  detectors  often  contain  a  multi¬ 
tude  of  different  electromagnetic  components.  Other 
than  the  use  of  high  and  low  dielectric  materials  men¬ 
tioned  above,  the  materials  can  be  isotropic  and  biax¬ 
ial  electromagnetically,  absorbing  and  non-absorbing, 
metallic  and  non-metallic,  and  form  homogenous  and 
periodic  layers  within  the  structure.  Modeling  of  such 
devices  therefore  requires  one  to  be  able  to  handle  suc¬ 
cessfully  all  these  material  components. 

The  electromagnetic  properties  of  such  a  complex 
multilayered  structure  depend,  sometimes  sensitively, 
on  a  large  number  of  material  parameters  such  as  the 
dielectric  constants  of  the  various  components  along 
different  directions,  and  the  widths  and  thicknesses  of 
individual  components  within  each  layer.  In  detector 
design,  this  set  of  parameters  must  be  optimized  under 
a  variety  of  constraints  to  achieve  a  certain  performance 
goal.  Modeling  of  photodetectors  therefore  requires  the 
use  of  methods  that  are  not  only  reliable,  accurate  and 
stable,  but  are  highly  efficient  and  have  only  modest 
memory  requirements. 

For  this  purpose,  we  have  developed(Jiang  et  al. 
2001)  a  modeling  method  based  on  the  determination 
of  the  modal  fields  within  each  layer  of  the  detector 
and  then  using  them  to  represent  the  actual  fields,  and 
the  use  of  the  transmission-line  technique  to  solve  the 
resulting  problem  associated  with  the  electromagnetic 
boundary  conditions.  Our  method,  referred  to  as  the 
modal  transmission-line  (MTL)  method,  is  fully  capa¬ 
ble  of  handling  the  stringent  requirements  mentioned 
above  in  connection  with  photodetector  modeling. 

In  the  MTL  method,  the  modal  fields  within  each 
of  the  periodic  layer  are  computed  from  the  eigenvalues 
and  eigenvectors  of  a  secular  equation  associated  with 
Maxwell’s  equations.  Our  work  here  is  to  describe  a 
method  to  transform  such  an  equation  into  a  form  that 
can  be  solved  very  efficiently  in  the  case  of  a  normally 
incident  wave.  We  then  apply  our  method  to  illustrate 
how  it  can  be  used  in  modeling  the  C-QWIP  mentioned 
above. 

Throughout  this  paper,  we  assume  that  the  inci¬ 
dent  field  propagates  perpendicularly  to  the  structure, 
and  that  the  entire  structure  can  be  decomposed  into 
a  multilayered  one,  where  a  common  unit  cell  can  be 
identified  for  all  the  periodic  layers.  This  means  that 
the  periodicities  of  each  of  the  periodic  layers  must  be 
commensurate  with  each  other.  We  further  assume 


that  within  each  periodic  layer,  the  unit  cell  exhibits 
reflection  symmetry  about  a  constant-a;  plane  passing 
through  its  center. 

Symmetry  is  a  familiar  concept  in  the  study  of 
electron  waves  in  crystalline  materials  and  photonic 
bandgap  crystals.  Symmetries  in  electromagnetic 
waves  are  more  interesting  because  other  than  the  sym¬ 
metries  associated  with  the  structure  itself,  there  are 
internal  symmetries  due  to  Maxwell’s  equations  them¬ 
selves. 

In  the  case  of  lamellar  structures  with  a  simple  unit 
cell,  the  symmetry  turns  out  to  be  rather  simple,  as  we 
will  discuss  in  this  work.  First  we  consider  the  case 
where  the  incident  wave  is  symmetric  just  like  the  unit 
cell.  An  obvious  example  is  the  very  common  case  of  a 
normally  incident  plane  wave.  Owing  to  the  symmetry 
of  the  unit  cell,  the  modal  fields  in  the  structure  can  be 
classified  as  either  symmetric  or  anti-symmetric.  Since 
the  incident  fields  are  symmetric,  they  do  not  couple 
to  the  anti-symmetric  modes  at  all.  By  explicitly  elim¬ 
inating  such  irrelevant  modes  in  modeling  symmetric 
structures  within  our  modal-field  representation  of  the 
scattering  problem,  we  find  that  a  numerical  computa¬ 
tional  speedup  of  up  to  8  times  can  be  achieved.  On 
the  other  hand,  if  the  incident  wave  is  anti-symmetric 
with  respect  to  the  symmetry  plane,  then  only  the  anti¬ 
symmetric  modal  fields  are  coupled.  An  up  to  8  times 
speedup  can  also  be  achieved. 

Symmetrized  modal  fields  can  be  employed  as  long 
as  the  incident  wave  falls  normally  onto  the  structure. 
Such  an  approach  offers  a  significant  numerical  advan¬ 
tage  even  when  the  incident  wave  does  not  possess  a 
definite  symmetry.  The  reason  is  that  a  general  nor¬ 
mally  incident  wave  can  always  be  decomposed  into  a 
symmetric  part  and  an  anti-symmetric  part.  In  addi¬ 
tion,  linear  combinations  of  the  Bloch  plane  waves  (i.e. 
plane  waves  that  obey  Bloch’s  theorem)  that  exhibit 
symmetric  and  anti-symmetric  properties  can  also  be 
made.  Owing  to  the  linearity  of  the  scattering  problem, 
these  two  different  types  of  waves  can  be  dealt  with  sep¬ 
arately  using  the  present  approach.  A  speedup  of  up 
to  4  times  can  still  be  achieved  over  the  conventional 
approach  which  ignores  symmetry  considerations. 

After  deriving  the  resulting  secular  equations  for 
the  symmetric  and  anti-symmetric  modes  for  the  TE 
and  TM  cases,  we  apply  our  symmetrized  modal  field 
method  to  compute  the  scattering  properties  of  the 
above  C-QWIP  structure.  We  also  compare  and  con¬ 
trast  the  results  with  those  obtained  from  conventional 
modal-field  method.  First  of  all,  we  find  that  the 
numerical  results  for  the  eigenvalues  and  eigenvectors 
from  two  methods  are  practically  identical  within  ma¬ 
chine  accuracy.  The  agreement  in  the  final  result  for 
the  absorption  is  better  than  1  x  10— 11 ,  but  the  results 
are  obtained  much  quicker  using  the  new  method.  The 
speedup  factor  indeed  approaches  closely  the  theoreti- 
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Therefore  one  has 


cal  expected  value  of  8  as  the  number  of  harmonics  in¬ 
cluded  in  the  calculation  is  increased.  Furthermore  the 
present  method  reduces  the  required  amount  of  mem¬ 
ory  by  about  one  half. 

2  SYMMETRIZED  MODAL-FIELD 
METHOD 

We  are  interested  in  a  biaxial  lamellar  structure 
that  is  periodic  along  the  x-direction.  If  a  is  the  size  of 
the  primitive  cell,  then  the  reciprocal  lattice  vectors 
along  the  x-direction  have  magnitudes  K  =  2Trn /a, 
where  n  =  0,±1, ....  The  basic  structure  of  inter¬ 
est  here  is  exactly  the  same  as  that  studied  in  the 
past,(Tamir  and  Zhang  1996,  Jiang,  Tamir  and  Zhang 
2001)  except  for  the  additional  requirement  that  the 
structure  must  be  symmetric  about  the  center  of  each 
primitive  cell.  Thus  if  the  origin  in  x  is  chosen  to  be 
at  the  center  of  one  of  the  primitive  cells,  the  dielectric 
constant  satisfies 

eQ(— x)  =  e“(+x),  for  a  =  x,y,z.  (1) 

We  therefore  consider  some  general  results  for  Fourier 
expansions  and  series  that  will  be  needed  later.  First 
any  function  /(x)  having  the  same  periodicity  as  the 
lattice  can  be  expanded  as  a  Fourier  series: 

f(x)  =  Y  f*  exp (iKx)  (2) 

X 

where  the  sum  over  K  implies  a  sum  over  the  integers 
that  specify  K ,  and  the  Fourier  coefficients  Jk  are  given 

by 

W2  dx 

fx  =  —  f(x)  exp(-iKx).  (3) 

J  —  a/2  a 

Alternatively  /(x)  can  be  expressed  in  the  form  of  a 
trigonometric  Fourier  series 

f(x)  =  fo+Y  [fK  cos(Kx )  +  f-x  sin (Kx)  .  (4) 

K>0 

It  is  easy  to  see  that  these  two  sets  of  Fourier  coefficients 
are  related  by 

fo  =  fo,  fx  =  fx  +  f-x,  f-x  =  i(fx  ~  f-x)-  (5) 

The  only  restriction  on  /(x)  thus  far  is  the  periodic 
property. 

Next  if  /(x)  is  a  symmetric  function  under  reflec¬ 
tion  about  the  x  =  0  plane,  i.e.  f{—x)  =  /(x),  so  that 
/(x)  is  an  even  function  of  x.  In  that  case  one  can  eas¬ 
ily  see  from  Eq.  (3)  that  /_ k  =  fx-  Using  this  relation 
in  Eq.  (2)  one  finds  that 

f(x)  =  fo  +  2  Y]  fx  cos(Kx).  (6) 

x>o 


fo  =  fo,  fx  =  2fx,  f-x  =  0,  (7) 

which  can  also  be  deduced  from  Eq.  (5).  For  our  calcu¬ 
lation  later,  we  find  that  it  sometimes  more  convenient 
to  rewrite  the  Fourier  series  expansion  of  a  symmetric 
function  /(x)  in  the  form 

f(x)  =  fo+Y  ^x  cos(Kx)  +  Y  f~x  cos  (Kx) 
x>o  x<o 

=  Y  fK  C0S(Kx )»  (8) 

x 

where  the  sum  over  K  is  unrestricted.  We  will  next 
derive  expressions  for  the  modal  fields.  We  find  that  it 
is  more  convenient  to  use  Eq.  (4)  to  express  the  modal 
fields,  but  Eq.  (8)  is  better  for  the  dielectric  constants. 

3  MODAL  FIELDS 

We  consider  here  the  modal  fields  inside  a  sym¬ 
metric  periodic  layer  located  between  z  =  0  and  t. 
The  incident  wave  is  assumed  to  be  propagating  along 
the  z-direction,  i.e.  normal  to  the  layer.  We  want  to 
derive  a  secular  equation  that  determines  the  modal 
fields  within  the  periodic  layer.  The  general  fields  are 
then  expressed  in  terms  of  these  modal  fields.  (Tamir 
and  Zhang  1996,  Jiang  et  al.  2001) 

3.1  TE  Mode 

In  the  case  of  TE  waves,  the  y-component  of  the 
electric  field  obeys  the  equation 

[k20ex(x)  +  dl]Ey  =  -d2zEy.  (9) 

In  previous  studies,  Ey  is  expressed  in  terms  of  plane 
waves  as 

Ey=Y  exp(iA'x)  E  Q‘X,m'Vm(z')  >  (10) 

X  m 

where 

Vm{z)  =  fm  exp (iKmz)  +  gm  exp (mm(i  -  z))  (11) 

are  the  equivalent  voltages,  m  labels  the  modes,  and 
the  transverse  component  of  the  wavevector  has  been 
set  to  zero.  It  is  known  that  the  coefficients  a,K,m  and 
Km  are  related  to  the  eigenvectors  and  eigenvalues  of 
the  secular  equation 

E  [(^,2  +  Km)dx,X'  —  &Q £x-x']ax',m  =  0i  (12) 

X' 

where  Sk,x'  is  the  Kronnecker  delta  function.  If  a  fixed 
integer  N  is  chosen  such  that  only  the  reciprocal  lat¬ 
tice  vectors  between  —2ttN/ci  and  2xN/a  are  retained 
in  our  calculation,  then  the  secular  equation  becomes 
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an  eigenvalue-eigenvector  problem  involving  square  ma¬ 
trices  of  size  (2 TV  +  1)  by  (2 N  +  1). 

It  is  important  to  note  that  a  symmetric  wave, 
such  as  a  plane  wave,  incident  onto  a  symmetric  struc¬ 
ture  can  only  couple  to  modes  that  are  also  symmetric. 
Therefore  instead  of  using  a  Fourier  series  expansion  of 
the  fields,  which  exhibit  no  definite  symmetry,  it  is  bet¬ 
ter  to  use  trigonometric  Fourier  series.  Then  we  only 
need  to  include  terms  that  are  symmetric,  i.e.  the  con¬ 
stant  term  and  the  cosine  terms.  Thus  in  Eq.  (9)  we 
let 

Ev  =  E  cos(Kx)  E  (*),  (13) 

K>  0  m 

and  use  the  cosine  expansion  of  the  form  given  in  Eq. 
(8)  for  the  dielectric  function 

ex(x)  =  £  cos (Kx).  (14) 

K 

Then  we  multiply  the  resulting  equation  by  cos(Kx), 
with  K  >  0,  and  integrate  within  a  unit  cell  to  obtain 
the  secular  equation 

Ck(K  T  ^m)^K,m  k0  ^  ]  Cft  k'Q'K'  ,m  =  9.  (15) 

K’>  0 


where  we  define 


Cic  = 

(2,  if  K  =  0 
(  1,  otherwise 

(16) 

and 

f+ 

tK,K' 

=  €K-K ’  +  f-K+K'  ■ 

(17) 

If  we  keep  only  values  of  K  between  0  and  2Nn/a, 
the  above  equation  can  be  solved  as  a  finite  order  ma¬ 
trix  eigenproblem  for  the  eigenvalues  and  eigenvectors, 
which  yield  and  a,K,m ,  respectively.  There  is  a 
total  of  N  +  1  of  these  symmetric  mode  eigenvalue- 
eigenvectors  pairs.  They  are  part  of  the  eigenvalue- 
eigenvectors  pairs  of  the  previous  secular  equation  [Eq. 
(12)].  We  have  checked  to  verify  that  this  is  indeed 
so.  Of  course  Eq.  (8)  also  contains  the  eigenvalue- 
eigenvectors  pairs  of  the  anti-symmetric  modes. 

A  secular  equation  for  the  anti-symmetric  modes 
can  be  derived  by  expanding  Ey  in  terms  of  only  the 
sine  function  rather  than  the  cosine  function: 

Ey  =  E  sin(A"x)  E  0’—K,mvm{z).  (13) 

K>  0  m 

Inserting  Eq.  (18)  in  Eq.  (9),  using  the  cosine  expansion 
of  the  form  given  in  Eq.  (8)  for  the  dielectric  function, 
multiplying  the  resulting  equation  by  sin(ATa;)  and  in¬ 
tegrating  x  over  a  unit  cell,  yields  the  secular  equation 
for  the  anti-symmetric  TE  modes: 

(A*  T  ^m)h— k0  ^  K'  =  H  (19) 

K’>  0 


In  most  practical  situations,  the  incident  field  is 
itself  symmetric  with  respect  to  the  symmetry  plane.  A 
very  common  example  is  that  of  a  plane  incident  wave. 
In  that  case,  the  anti-symmetric  modes  are  not  coupled 
at  all  to  the  incident  wave,  and  therefore  need  not  be 
computed.  To  obtain  the  modal  fields  within  each  of 
the  periodic  layers,  all  we  need  is  to  solve  the  secular 
equation  for  the  symmetric  modes  using  Eq.  (15).  By 
explicitly  dealing  with  the  relevant  symmetric  modes, 
we  now  only  need  to  solve  a,  N  +  1  by  N  +  1  matrix 
instead  of  a  2N  + 1  by  21V  + 1  matrix  for  its  eigenvalues 
and  eigenvectors.  Since  the  numerical  complexity  of 
matrix  eigenproblem  goes  as  the  third  power  of  the  size 
of  the  matrix,  we  see  that  the  computation  should  be 
about  8  times  faster,  since  N  is  typically  a  large  integer. 

On  the  other  hand,  in  the  rare  case  where  the  in¬ 
cident  field  is  anti-symmetric,  then  we  need  to  solve 
the  secular  equation  for  the  anti-symmetric  modes  us¬ 
ing  Eq.  (19).  We  now  only  need  to  solve  an  N  by 
N  matrix  instead  of  a  2 N  +  1  by  2N  +  1  matrix  for 
its  eigenvalues  and  eigenvectors.  Therefore  we  can  also 
achieve  a  speed-up  factor  of  about  8. 

3.2  TM  Mode 

In  the  case  of  TM  waves,  the  y-component  of  the 
magnetic  field  obeys  the  equation 

[k%ex( x )  +  ex(x)dxrjz(x)dx]Hy  =  - d\Hy ,  (20) 

where  i]z(x)  =  l/ez(x).  In  previous  studies,  Hy  is  ex¬ 
pressed  in  terms  of  plane  waves  as 

hv  =  Y,  exp(iKx)  E  bK,mim(z),  (21) 

K  m 

where 

im(z)  =  fm  exp (iKmz)  -  gm  exp (mm{t  -  z)),  (22) 

are  equivalent  currents.  The  coefficients  bK,m  and  Krn 
are  related  to  the  eigenvectors  and  eigenvalues  of  the 
secular  equation 

E  (^0 eK-K’  —  ^mSK,K')bK\m  (23) 

K' 

—  E  E  Ek-K'E  'VK'-K”E,,bK”,m  =  0. 

K'  K” 

Equation  (23)  expresses  a  matrix  eigenproblem  for  a 
2 N  +  1  by  2N  +  1  matrix  if  reciprocal  lattice  vectors 
between  —  2nN/a  and  2nN/a  are  retained. 

It  is  known  that  (Li  1996)  faster  convergence  is 
achieved  if  e^_Ki  is  obtained  by  first  inverting  the 
matrix  containing  the  Fourier  coefficients  of  1  / f\_K, 
and  then  taking  the  K K'  element.  Similarly  gzK_Ki  is 
obtained  by  first  inverting  the  matrix  containing  the 
Fourier  coefficients  of  ezK_K,  and  then  taking  the  KK' 
element. 


4 


We  now  derive  the  secular  equation  for  symmetric 
TM  modes  in  the  trigonometric  Fourier  series  represen¬ 
tation  by  including  in  the  expansion  for  the  magnetic 
field  only  symmetric  terms 

hv=Y,  cos(Kx)  E  bK,m,im{z),  (24) 

K>  0  m 

and  use  the  cosine  expansion  of  the  form  given  in  Eq. 
(8)  for  the  dielectric  function  and  an  analogous  one  for 
r]z  ( x ) 

Vz{x)  =  Y  V k  cos (Kx).  (25) 

K 

We  use  Eq.  (24),  Eq.  (14),  and  Eq.  (25)  in  Eq.  (20), 
and  then  multiply  the  resulting  equation  by  cos  (A' a:), 
with  K  >  0,  and  integrate  within  a  unit  cell.  After 
some  algebra  we  obtain  the  following  secular  equation 
for  the  symmetric  TM  modes 

Y  (C K^m^K,K'  ~  (26) 

K'>  0 

+  ^  6K,K'K  ^k' ,K”K”bK” ,m  ~  0- 

K'> 0  K”> 0 

This  secular  equation  does  not  contain  the  unwanted 
anti-symmetric  modes  and  as  a  result  its  eigenvalues 
and  eigenvectors  can  be  computed  about  8  times  faster 
than  if  one  uses  Eq.  (23) . 

To  find  the  secular  equation  for  the  anti-symmetric 
modes,  we  use  the  sine  expansion  for  the  magnetic  field 

Hy=Y  sin  (Kx)  E  b—K,m^m  (z).  (27) 

K>  0  m 

Using  Eq.  (27),  Eq.  (14),  and  Eq.  (25)  in  Eq.  (20), 
and  multiplying  the  resulting  equation  by  sin(A^a:), 
with  I\  >  0,  and  then  integrating  within  a  unit  cell, 
we  obtain  the  following  secular  equation  for  the  anti¬ 
symmetric  TM  modes 

^rrif—K,m  k0  E  ^K,K'  b—K’  ,m  (28) 

K’>0 

+  Yj  E  eK,K'^-  =  0. 

K’>0  K”>0 

For  a  truncation  value  of  TV,  this  secular  equa¬ 
tion  become  an  TV  by  TV  matrix  eigenproblem  whose 
eigenvalues  and  eigenvectors  determine  the  modal  fields 
within  such  a  periodic  layer.  In  the  case  of  an  anti¬ 
symmetric  incident  wave,  we  can  use  this  secular  equa¬ 
tion  to  find  the  modal  fields  about  8  times  faster  than 
using  the  conventional  secular  equation. 

4  ALTERNATIVE  DERIVATION 

The  above  secular  equations  for  the  symmetric  TE 
and  TM  modes  can  also  be  derived  by  considering  the 


relations  between  the  coefficients  for  the  fields  in  the 
Fourier  series  representation  with  those  in  the  trigono¬ 
metric  Fourier  series  representation.  The  general  rela¬ 
tion  is  given  by  Eq.  (5).  The  derivation  of  the  secular 
equations  is  more  convenient  if  matrix  notation  is  used. 

For  that  purpose  we  define  the  following  column 
vectors 

'  h  i  r  /-i 

fK=  h  ,  f-K  =  ,  (29) 


'  A  i  r  l-i  ' 

h  =  h  ,  f-K  =  f~  2  •  (30) 


We  also  define  composite  vectors 


/o 

fo 

/= 

fh' 

,  /  = 

fh' 

f-K 

.  f-K  . 

and  the  transformation  matrix 


U  J.  J.  , 

L  o  a  —H 


(31) 


(32) 


where  1  is  an  identity  matrix.  The  relation  given  by 
Eq.  (5)  can  then  be  expressed  as 

/  =  Sf.  (33) 

The  inverse  relation  is  given  by 

f  =  S~1f,  (34) 

where 

'10  O' 

S”1  =  0  ±1  -§1  .  (35) 

.  0  \1  |1 

We  define  vectors  ex  and  €-k  analogous  to  Eq.  (31), 
and  a  matrix  £k,k  whose  ij  element  is  given  by  e k,  -  k.  ■ 


4.1  TE  Modes 

We  reconsider  first  the  TE  modes.  In  matrix  nota¬ 
tion  we  see  that  the  secular  equation  in  Eq.  (12)  in  the 
usual  plane  wave  basis  can  now  be  written  as 

(e  —  D)a  =  0,  (36) 

where  we  have  defined 


eo 

eT 

eK 

i 

tc 

h 

UJ 

*  (37) 

£k 

eK,K 

eK-K 

£-k 

e 

-K,K 

e-K-K 

'  K2 

0 

0  ' 

0 

r 

0 

(38) 

0 

0 

r 

T  =  K2  +  k21 , 


(39) 

where  K  is  a  diagonal  matrix  with  diagonal  elements 
K i,  K2, . . .,  and  the  composite  vector  a  is  defined  as  in 
Eq.  (31).  From  Eqs.  (34)  one  has  a  =  S-1a  and  so 
we  see  that  the  secular  equation  in  the  trigonometric 
Fourier  series  basis  is 

(e  -  D)a  =  0,  (40) 

where  for  any  2 N  +  1  by  2 TV  +  1  matrix  M ,  we  define 
the  matrix 


that  these  two  parts  can  be  treated  separately  using 
their  respective  secular  equations  within  the  periodic 
layer.  Consequently  we  expect  to  have  a  speed-up  fac¬ 
tor  of  up  to  4  for  high  truncation  value,  N. 

4.2  TM  Modes 

Next  we  reconsider  the  TM  case  using  this  matrix 
approach.  The  secular  equation  in  the  plane  wave  ba¬ 
sis  is  given  by  Eq.  (23),  which  can  now  be  written  in 
matrix  notation  as 


M  =  SMS-1, 


(41) 


(e  —  k21  -  eKf/K)b  =  0, 


(45) 


which  is  related  to 
One  finds  that 


e  = 


eo 

2eK 

0 


ca 

e+ 

fcA,A 

0 


M  by  a  similarity  transformation. 
D  D  and 


0 

0 

eK,K  . 


(42) 


where  we  have  defined 


CA,A 


—  cat, a-  ±  £k,-k- 


(43) 


Hence,  we  can  prove  that  the  secular  equation  for  TE 
waves  in  the  trigonometric  Fourier  series  basis  is  given 

by 


e0  -  a2  e-K  0 

2eK  e+  K  -  T  0 
0  0  eK,K  ~  r 


X 


X 


a0 

aK 

a~K 


=  0. 


(44) 


The  eigenproblem  is  seen  to  decouple  into  two  separate 
diagonal  blocks.  The  first  block,  consisting  of  the  first 
two  equations,  is  for  the  symmetric  modes.  It  is  easy 
to  see  that  they  are  equivalent  to  the  secular  equation 
that  we  derived  earlier  in  Eq.  (15). 

The  second  block,  consisting  of  the  last  equation, 
is  clearly  the  secular  equation  for  the  anti-symmetric 
modes.  This  equation  is  in  the  form  of  an  eigenproblem 
of  an  N  by  N  matrix.  If  the  incident  wave  is  symmet¬ 
ric,  these  anti-symmetric  modes  cannot  be  excited  and 
therefore  there  is  no  need  to  solve  this  secular  equation 
at  all.  However  if  the  incident  wave  is  anti-symmetric, 
then  we  need  to  solve  for  the  eigenvalues  and  eigen¬ 
vectors  of  this  secular  equation,  but  not  those  for  the 
symmetric  modes.  We  also  get  a  speed-up  factor  of  up 
to  8  times  compared  with  the  conventional  approach. 

It  is  important  to  note  that  our  approach  of  using 
symmetrized  modal  functions  is  applicable  even  when 
the  incident  wave  has  no  definite  symmetry.  The  reason 
is  that  a  general  incident  wave  can  always  be  written 
as  the  sum  of  a  symmetric  part  and  an  anti-symmetric 
part.  Linearity  of  the  electromagnetic  problem  implies 


where  1  is  an  identity  matrix,  and 


Vo 

T 

Vk 

T 

V-K 

77a 

77a,  A 

7?  A, -A 

77-a 

»7— AT.A 

7?- A, -A 

(46) 


This  equation  can  then  be  transformed  into  the  trigo¬ 
nometric  Fourier  series  basis  by  multiplying  from  the 
left  by  S.  The  result  is 

—  k21  —  eKfjK^j  6  =  0,  (47) 

where 


1  =  1, 


and 


V  = 


VO  V  A  0 

2t?a  Vk,  k  0 

0  0  r]KK 


with 


Vk,k  =  Vk,k 


±  77a- A- 


(48) 


(49) 


(50) 


We  can  easily  compute  K  and  obtain  from  Eq.  (47) 
the  following  secular  equation  for  the  TM  modes  in  the 
trigonometric  Fourier  series  basis 


/  r  e0  0 

2ca  e+A^-  0 

VL  0  0  eK.K$r 


X 


X 


0, 


where 


«*  =  (!  -Kt,±'KK)' 

and 


I  = 


10  0 
0  10 
0  0  1 


(51) 


(52) 


(53) 
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Figure  2:  The  absorption  at  a  fixed  wavelength  (near 
the  maximum)  is  plotted  as  a  function  of  the  truncation 
factor  N. 


Figure  3:  The  absorption  spectrum  computed  using  a 
fixed  value  of  the  truncation  factor  N. 


is  the  identity  matrix.  The  matrix  in  Eq.  (51)  is 
seen  to  separate  into  two  diagonal  blocks.  The  first 
block,  which  consists  of  two  parts,  clearly  determines 
the  modes  of  the  symmetric  TM  modes.  The  result 
can  be  seen  to  be  equivalent  to  the  secular  equation 
that  we  derive  earlier  in  Eq.  (26).  The  second  block 
on  the  other  hand  determines  the  modes  of  the  anti¬ 
symmetric  TM  waves.  These  modes  are  not  needed  for 
symmetric  normally  incident  waves. 

Our  alternate  derivation  of  the  secular  equation 
clearly  shows  that  both  the  symmetric  and  the  anti¬ 
symmetric  modes  are  contained  in  the  conventional 
method  that  ignores  symmetry  considerations.  The  use 
of  symmetrized  modal  field  here  enables  one  to  ignore 
modes  that  do  not  have  the  right  symmetry  and  are 
therefore  irrelevant  to  the  problem,  and  hence  leads  to 
in  a  more  efficient  numerical  method. 

5  RESULTS  AND  DISCUSSION 

We  will  show  the  results  of  numerical  computations 
of  the  scattering  properties  of  a  symmetric  structure  at 
normal  incidence  using  the  present  method  and  com¬ 
pare  them  with  those  obtained  using  the  conventional 
method.  Our  purpose  here  is  to  illustrate  the  superior¬ 
ity  of  the  present  method  by  considering  as  an  example 
a  practical  C-QWIP  structure,  (Choi  1997,  Mao  et  al. 
2002,  Choi  et  al.  2003,  2004a, b,c)  as  shown  in  Fig.  1. 

First  of  all,  we  verify  that  the  eigenvalues  and 
eigenvectors  obtained  from  the  present  method  and  the 
conventional  method  agree.  Indeed  we  find  that  they 
are  practically  identical  within  machine  accuracy.  One 
of  the  most  important  quantity  of  interest  in  modeling 
QWIPs  is  the  relative  absorption,  P;i\)S  at  any  given 
wavelength.  We  find  that  using  the  resulting  eigenval¬ 
ues  and  eigenvectors  to  compute  Pakg,  the  agreement 


is  better  than  1  x  10-11. 

Typically  in  our  modeling  work,  the  absorption 
spectrum  is  first  computed  using  a  moderate  value  of 
the  truncation  factor,  N.  To  test  for  convergence  of 
the  results,  we  then  fix  the  wavelength  near  one  of  the 
absorption  peaks  and  repeat  the  computation  of  the 
spectrum  using  either  methods  for  larger  values  of  N. 
As  a  result  we  obtain  a  convergence  curve  as  shown 
in  Fig.  2.  This  curve  provides  us  a  good  idea  of  the 
maximum  error  in  the  computed  result  at  any  given  N. 
Using  a  value  of  N  large  enough  to  ensure  that  our  re¬ 
sults  are  accurate  within  a  certain  acceptable  tolerance, 
we  then  recompute  the  spectrum.  The  result  is  shown 
in  Fig.  3. 

The  above  results  were  computed  using  the  pre¬ 
sent  symmetrized  modal  field  method  as  well  as  the 
conventional  method.  Since  they  agree  to  better  than 
1  x  1CU11  for  any  wavelength  and  truncation,  the  re¬ 
sults  from  the  above  graphs  will  be  practically  identi¬ 
cal.  However  the  present  method  reduces  the  memory 
requirement  to  about  one  half  that  of  the  conventional 
method.  It  also  runs  much  faster,  as  can  be  seen  from 
Fig.  4,  where  the  computation  times  for  obtaining  the 
absorption  spectrum  using  these  two  methods  are  plot¬ 
ted  as  a  function  of  N.  We  define  a  speed-up  factor 
as  the  ratio  of  the  computation  time  using  the  present 
method  to  that  of  the  conventional  method.  The  speed¬ 
up  factor  is  shown  in  Fig.  5  as  a  function  of  N.  We  see 
that  the  speed-up  factor  increases  with  N  and  indeed 
gradually  approaches  the  expected  value  of  8. 

For  the  present  somewhat  simple  C-QWIP,  rather 
accurate  absorption  spectrum  can  be  obtained  in  about 
11.4  seconds  for  a  given  structure  using  a  truncation 
value  of  N  =  15.  However  in  optimizing  the  perfor¬ 
mance  of  the  detector,  even  if  the  material  parame¬ 
ters  are  kept  fixed,  one  still  needs  to  explore  numerous 
structures  with  layers  having  different  widths,  spacings 
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Figure  4:  The  total  computation  time  for  obtaining 
an  absorption  spectrum  is  plotted  as  a  function  of  the 
truncation  factor  N. 
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Figure  5:  The  speed-up  factor  as  a  function  of  the  trun¬ 
cation  factor  N. 


and  thicknesses.  For  example  if  we  want  to  vary  only  4 
out  of  the  total  of  6  dimensional  parameters,  and  if  each 
parameter  can  assume  only  10  different  values,  then 
an  exhaustive  exploration  of  detector  performance  will 
take  close  to  32  hours  using  the  conventional  method. 
The  speed-up  factor  for  the  new  method  at  this  trunca¬ 
tion  value  is  about  3.2.  That  means  that  the  modeling 
can  now  be  done  in  only  about  10  hours.  In  general, 
the  saving  will  be  much  greater  for  more  complicated 
structures (Gunapala  et  al.  2000)  and  for  higher  accu¬ 
racies. 
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